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ABSTRACT 

Context. High-mass stars play a prominent role in Galactic evolution, but their formation mechanism is still poorly understood. This 
lack of knowledge reflects the observational limitations of present instruments, whose angular resolution (at the typical distances 
of massive protostars) precludes probing circumstellar gas on scales of 1-100 AU, relevant for a detailed investigation of accretion 
structures and launch/collimation mechanims of outflows in high-mass star formation. 

Aims. This work presents a study of the milliarcsecond structure of the 6.7 GHz methanol masers at high-velocity resolution 
(0.09 km s" 1 ) in four high-mass star-forming regions: G16.59-0.05, G23.01-0.41, IRAS 20126+4104, and AFGL 5142. 
Methods. We studied these sources by means of multi-epoch VLBI observations in the 22 GHz water and 6.7 GHz methanol masers, 
to determine the 3-D gas kinematics within a few thousand AU from the (proto)star. Our results demonstrate the ability of maser 
emission to trace kinematic structures close to the (proto)star, revealing the presence of fast wide-angle and/or collimated outflows 
(traced by the HtO masers), and of rotation and infall (indicated by the CH3OH masers). The present work exploits the 6.7 GHz maser 
data collected so far to investigate the milliarcsecond structure of this maser emission at high-velocity resolution. 
Results. Most of the detected 6.7 GHz maser features present an ordered (linear, or arc-like) distribution of maser spots on the plane 
of the sky, together with a regular variation in the spot LSR velocity (Vlsr) with position. Typical values for the amplitude of the Vlsr 
gradients (defined in terms of the derivative of the spot Vlsr with position) are found to be 0.1-0.2 km s mas~'. In each of the four 
target sources, the orientation and the amplitude of most of the feature Vlsr gradients remain remarkably stable in time, on timescales 
of (at least) several years. We also find that the data are consistent with having the V LSR gradients and proper motion vectors in the 
same direction on the sky, considered the measurement uncertainties. In three (G16.59-0.05, G23. 01-0.41, and IRAS 20126+4104) 
of the four sources under examination, feature gradients with the best determined (sky-projected) orientation divide into two groups 
directed approximately perpendicular to each other. 

Conclusions. The time persistency, the ordered angular and spatial distribution, and the orientation generally similar to the proper 
motions, altogether suggest a kinematical interpretation for the origin of the 6.7 GHz maser Vlsr gradients. This work shows that 
the organized motions (outflow, infall, and rotation) revealed by the (22 GHz water and 6.7 GHz methanol) masers on large scales 
(~ 100-1000 AU) also persist to very small (~10 AU) scales. In this context, the present study demonstrates the potentiality of the 
mas-scale 6.7 GHz maser gradients as a unique tool for investigating the gas kinematics on the smallest accessible scales in proximity 
to massive (proto)stars. 

Key words. Masers - Techniques: high angular resolution - Techniques: spectroscopic - ISM: kinematics and dynamics - ISM: struc- 
ture - Stars: formation 

1. Introduction iMoscadelli et alj [201 lb . and they have also been successfully 

used t o infer the propertie s of the magnetic field inside molecular 
Intense maser transitions of the OH (at 1.6 and 6.0 GHz), cores (|Fi sh & Reicj | 2006 | ) . 

CH3OH (at 6.7 and 12 GHz), H2O (at 22 GHz), and SiO VLBI maser observations are useful for testing the present 
(at 43 GHz) molecules, are observed towards massive star- theories of massive star formation. To overcome the radiation 
forming regions. The elevate brightness temperature (>10 7 K) pressure exerted by the already ignited star and still enable mass 
of the maser emissions permits us to observe them with the accretion, most models point to the essential role that accre- 
Very Long Baseline Interferometry (VLBI) technique, achiev- tion disks and jets should play. The former focus the ram pres- 
ing both very high angular (1-10 mas) and velocity (0. 1 km s _1 ) sure of the accreting matter across the disk plane, and the latter 
resolution. Maser VLBI observations are the unique mean by channel stellar photons along the jet axis and lower the radia- 
which one can explore the gas kinematics close (within tens tion pressure across the equatorial plane. However, so far, in- 
or hundreds of AU) to the forming high-mass (proto)star (e.g. , terferometric observations of thermal continuum and line emis- 
iGoddi etalJl2005t Sanna et a D 1201 Pair* iMatthews et al. 2010; sions towards massive star-forming regions have only identi- 
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fied a handful of disk candidates i n association with high-mass 
stars (see e.g JCesaroni et al.l2006l) . In particular, no circumstel- 
lar disks have been detected around early O-type stars, where 
only huge, massive, rotating structures are seen, whose lifetimes 
appe ar to be less than the corresponding rotation periods (see 
e.g. jBeltran et al.|[201 ll) . One reason that could explain the diffi- 
culty detecting disks in massive star-forming regions (with typi- 
cal distances of a few kpc) is that the angular resolution (> 1") 
achieved by present millimeter interferometers could not suffice 
to resolve disks if their size, as models predict, is comparable 
to or less than a few thousand AU. On the contrary, supposing 
one can observe a masing transition from the disk gas, the angu- 
lar resolution achievable with maser VLBI observations is high 
enough to accurately measure the rotation curve of the accretion 
disk. 

For a decade we have been studying a sample of ten can- 
didate high-mass (proto)stars by means of VLBI observations 
of H 2 22 GHz, CH 3 OH 6.7 GHz, and OH 1.6 GHz masers 
associated with the (proto)stellar environment. The H2O and 
CH3OH masers are observed at several epochs to derive the 
maser proper motions, which, knowing the distance to the 
maser source, can be translated into the sky-projected com- 
ponents of the maser velocity. Combining this information 
with the line-of-sight velocity components, which are derived 
from the maser LSR velocity (Vlsr) and the knowledge of 
the (proto)star systemic Vlsr, a full 3-D picture of the mo- 
tion of the masing gas can be obtained. S o far, data for four 
sources have been analyzed: G16 .59-0.05 (Sanna et al.l l2010a. 



hereafter SMC1), G23.01-0.41 dSanna et al l 1201 Obi hereafter 



SMC2), IRAS 20126+410 4 dMoscadelli et alJ 1201 ll hereafter 
MCR), and AFGL 5142 (Godd i et all 120071 I20TTI hereafter 
GMS1 and GMS2, respectively). These results demonstrate the 
ability of maser emission to trace kinematic structures close to 
the (proto)star, revealing fast wide-angle and/or collimated out- 
flows (traced by the H2O masers) and rotation and infall (indi- 
cated by the CH3OH masers). 

The present work exploits the 6.7 GHz maser data collected 
so far to investigate the milliarcsecond structure of this maser 
emission at high-velocity resolution. Since the velocity resolu- 
tion of the 6.7 GHz maser dataset (0.09 km s _1 ) is high enough 
to resolve the maser emission linewidth (with typical FWHM 
of 0.3 km s ), by mapping each velocity channel across the 
maser linewidth, one can investigate how the spatial structure 
of the emission changes with the velocity. A similar study has 
been previously per formed for the intense 1 .6 GHz dFish et al.l 
2006) and 6.0 GHz dFish & Siouwermanl2007l) OH and 12 GHz 
CH3OH dMoscadelli et al.l l2003h masers in the UC Hn region 
W3(OH). These observations have revealed that it is quite com- 
mon to observe a regular variation (along a line or an arc) in the 
maser peak position with the Vlsr- In W3(OH), 1 .6 and 6.0 GHz 
OH masers and 12 GHz CH3OH masers are found to have simi- 
lar amplitudes of the Vlsr gradients (calculated by dividing the 
maser Vlsr linewidth by the path length over which the maser 
peak position shifts on the plane of the sky), varying in the range 
0.01-1 km s" 1 ALT 1 . 

This paper extends the study of the milliarcsecond struc- 
ture of the maser emission at high-velocity resolution to 
the 6.7 GHz CH3OH masers, presenting results in four dis- 
tinct massive star-forming regions: G16.59-0.05, G23.01-0.41, 
IRAS 20126+4104, and AFGL 5142. Section|2]summarizes the 
main observational parameters of the EVN 6.7 GHz observations 
towards these sources. Section [3] describes the basic character- 
istics of the milliarcsecond, velocity structure of the 6.7 GHz 
masers, making use of our multi-epoch observations to also 



examine the time variation. Section [4] compares the directions 
on the sky plane of the maser Vlsr gradients and proper mo- 
tion vectors. Section|5]compares the results obtained among the 
four sources and discusses a possible interpretation of the time- 
persistent velocity gradients inside most of the 6.7 GHz maser 
cloudlets. 



2. Summary of CH 3 OH 6.7 GHz maser EVN 
observations 

Using the European VLBI Network (EVNfl we observed the 
sources G16.59-0.05, G23.01-0.41, IRAS 20126+4104, and 
AFGL 5142 in the 5i-6 A + CH 3 OH transition (rest frequency 
6.668519 GHz), at three different epochs, over the years 2004- 
2009. Table Q] lists the EVN observing epochs for each of the 
four sources, and gives the bibliographic reference of the articles 
where the results of these observations have been published. We 
refer the reader to these articles both for a full description of the 
EVN observing setup and for a discussion of the 3-D gas kine- 
matics close to the massive (proto)star as traced by the H2O 
22 GHz and CH3OH 6.7 GHz masers. 

To determine the 6.7 GHz maser absolute positions and 
velocities, the EVN observations were performed in phase- 
referencing mode, by fast switching between the maser source 
and one strong, nearby calibrator. For each of the four maser tar- 
gets, the selected phase-reference calibrator (belonging to the list 
of sources defining the International Celestial Reference Frame, 
ICRF) has an angular separation on the sky from the maser less 
than 3°. In each of the four targets, the absolute position of the 
6.7 GHz masers is derived with an accuracy of a few milliarc- 
seconds. 

The antennae involved in the observations were Medicina, 
Cambridge, Jodrell, Onsala, Effelsberg, Noto, Westerbork, 
Torun, Darnhall, and Hartebeesthoek. Because of technical prob- 
lems, the Darnhall telescope replaced the Jodrell antenna in 
the first epoch (November 2004) of the IRAS 20126+4104 and 
AFGL 5142 EVN observations. The Hartebeesthoek antenna 
took part in the observations only for the first two epochs of the 
low-DEC sources G16.59-0.05 and G23.01-0.41. In the third 
epoch, since the maser emission would be heavily resolved on 
the longest baselines involving the Hartebeesthoek antenna, it 
was replaced with the Onsala telescope. All the EVN runs had a 
similar observing setup, providing the same velocity resolution 
of 0.09 km s _1 , and a similar sensitivity, corresponding to a ther- 
mal rms noise level on the channel maps of ~ 10 mJy beam -1 . 

As indicated in Table [T] the EVN observations achieved 
comparable angular resolutions towards the four target sources, 
with some larger beams for the two sources at low DEC, because 
of the degraded Mv-coverage. Isolated 6.7 GHz maser emission 
centers on a single-channel map (named "spots") are generally 
unresolved or slightly resolved with the EVN. By fitting a two- 
dimensional elliptical Gaussian to the spot intensity distribution, 
the accuracy in the relative positions of maser spots, 56, can be 
calculated through the formula 56 = 0.5 FWHM/SNR, where 
FWHM is the fitted spot size (about the FWHM beam size for a 
compact source) and SNR is the ratio of the spot intensity with 
the image rms noise. Since most of the detected spots have inten- 
sities >1 Jy and SNR>100, their relative positions are known 
with an accuracy better than 0.1 mas. 



1 The European VLBI Network is a joint facility of European, 
Chinese, South African, and other radio astronomy institutes founded 
by their national research councils. 
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3. Results 

In the following we use the term "spot" to point to a compact 
emission center on a single-channel map and the term "feature" 
to refer to a collection of spots emitting in contiguous channels 
at approximately the same position on the sky (within the beam 
FWHM). In our view, a maser feature corresponds to a distinct, 
masing cloud of gas, whose spatial and velocity structure is the 
subject of our analysis. 

3.1. 6.7 GHz maser internal V L sr gradients 

Figure [Upresents the spatial distribution and the change in Vlsr 
with position of the spots belonging to two 6.7 GHz maser fea- 
tures, selected to be representative of the ordered Vlsr gradi- 
ents commonly found inside the 6.7 GHz masers. The exam- 
ples presented in Fig.Q]are two intense, persistent features from 
the source G23.01-0.41. At each of the three observing epochs, 
spots from a given feature distribute in space close to a line and 
show a regular variation in Vlsr versus the position measured 
along the feature major axis (see also Fig. 6 in SMC1). In the fol- 
lowing discussion, we use r s and r v to indicate the correlation 
coefficients of the linear fit to the spot positions on the sky plane 
and to the Vlsr variation with position (measured along the ma- 
jor axis of the spot distribution), respectively. For the maser fea- 
ture in Fig. [Tk, the linear fits of both the sky -projected spot dis- 
tribution and the variation in Vlsr with positions always have a 
correlation coefficient equal to 1. For the feature in Fig.QJ, the 
spatial and velocity distribution of the spots is less ordered, with 
correlation coefficients of the linear fits in the range 0.4-1.0. 

The linear fits presented in Fig. [1] were performed for each 
feature (of the four maser sources) with three or more spots, 
thereby deriving the following quantities: the P.A. of the feature 
(major axis) orientation on the sky plane, P s , with the degree of 
linear correlation of the spot distribution on the sky plane mea- 
sured by r s ; and the amplitude of the Vlsr gradient with posi- 
tion, F (in units of [km s mas -1 ]), with the degree of linear 
correlation of Vlsr with position (measured along the feature 
major axis) given by r v . With the adopted convention, that posi- 
tion's offsets along the feature major axis are taken as positive if 
increasing to east, Vlsr increasing to east (west) results into pos- 
itive (negative) velocity gradients. Then, the P.A. of the feature 
orientation, P s , varies in the range 0°-180° for positive gradi- 
ents, and in the range 180°-360° for negative gradients. The 
measurement errors of both the spot positions (< 0.1 mas) and 
Vlsr (^ 0.09 km s _1 ) are small enough with respect to the typi- 
cal extension of the feature spot distribution (from a few mas to 
10 mas) and feature emission width (0.5-1 km s _1 ) not to affect 
the precision with which the orientation and the amplitude of the 
feature gradient can be determined. 

The examples shown in Fig. Q] are not rare cases among the 
6.7 GHz maser features. For each of the four target sources, 
Fig. |2] presents the histograms of the correlation coefficients r s 
and r v of the linear fits performed on all the features with three 
or more spots. A large fraction (60%-80%) of the correlation co- 
efficients showing values higher than 0.6, indicates that, in each 
source, a vast majority of maser features has a spatial structure 
primarily elongated in one direction, and also presents (at least) 
a hint of a regular variation in Vlsr along this direction. Besides, 
more than half of the features present a well-defined linear struc- 
ture with quite a regular change in Vlsr with position, since 
more than half of the derived values of r s and r v are higher 
than 0.8. We note that, in each source, the distribution of r v is 
more peaked at higher values than that of r s . That means that 



there are features showing a good correlation of Vlsr with po- 
sition measured along the major axis of elongation, even if their 
structure is only marginally elongated. Figure QJ) illustrates this 
case for a representative feature. 

Figure [3] presents the histograms of the amplitude of the 
Vlsr gradients for the observed maser sources. To produce the 
histograms and calculate the average and standard deviation 
of the distributions, we selected only those features that pre- 
sented at least some evidence of being elongated (r s > 0.6) 
and that showed a hint of a regular Vlsr variation with position 
(r v > 0.6). The average amplitude of the Vlsr gradients for the 
sources G23.01-0.41 and G16.59-0.05 («0.2 km s _1 mas" 1 ) is 
about twice the value derived for the sources IRAS 20126+4104 
and AFGL 5142 (^0.1 kms" 1 mas" 1 ). 

The knowledge of an accurate value of the source distance 
would allow us to convert the amplitude of Vlsr gradients 
from units of [km s _1 mas -1 ] to units of [km s _1 AU~ ! ] 
and to compare the gradient amplitudes of different sources. 
Recently, the distance to two of the sources has been pre- 
cisely determined by measuring the parallax of the associated, 
strong methanol, and water masers: IRAS 20126+4104 is at 
a distance of 1.64+0.05 kpc (MCR ), and G23.01-0.41 is at 
4.59+0.4 kpc (Brunthal e'r et alj|2009b . After converting to units 
of [km s AU -1 ], one finds that the distribution of gradient am- 
plitudes in IRAS 20126+4104 (0.07+0.04 km s" 1 AIT 1 ) over- 
laps with the one in G23.01-0.41 (0.05+0.05 km s _I AIT 1 ), 
although the (proto)stars associated to the CH3OH masers in 
IRAS 20126+4104 and G23.0 1-0.41 have very different lumi- 
nosities and masses (MCR,SMC2). That lets us argue that the 
Vlsr gradients internal to the 6.7 GHz masers might have a 
common origin in the two sources, resulting in comparable val- 
ues of gradient amplitudes. Also for the sources G16.59-0.05 
and AFGL 5142 (using the more uncertain distance of 4.4 kpc 
and 1.8 kpc, respectively), the values of the gradient ampli- 
tude, in units of [km s" 1 AU -1 ], are consistent with those of 
IRAS 20126+4104 and G23.0 1-0.41. 

3.2. Time persistency of V L sr gradients 

As defined in Sect. 13.11 Vlsr gradients are vector quantities, ori- 
ented on the sky plane at RA. = P s and with amplitude T. 
For the two features presented in Fig.Q] the direction on the sky 
plane and the amplitude of the Vlsr gradient appear to be fairly 
constant over the three observing epochs. To investigate the vari- 
ation in Vlsr gradients with time, we selected features persistent 
over two or three epochs and calculated the standard deviation 
of the direction (A P s ) and amplitude (A F) of the feature Vlsr 
gradient at different epochs. The (time-average) mean gradient is 
calculated by taking the mean values of the gradient components 
projected along the east and the north directions. For a persistent 
feature, the change in time of the gradient amplitude is measured 
by the "fractional time variation" (AT / F), defined as the ratio 
of the standard deviation, A F, and the mean of the values of the 
feature gradient amplitude, F, at different epochs. 

Figure |4]presents the histograms of A P s and AT/T for 
each of the four maser sources. In each source, on a timescale 
of several years, the amplitude of the Vlsr gradient for most 
of the persistent features changes less than 20-40%, and the 
variation in the gradient orientation is generally less than 20°. 
Measurement errors can contribute to enlarge the calculated 
spreads in time, so that the effective time variation could be even 
smaller than what is shown in Fig. [4] Therefore, we can conclude 
that most of the 6.7 GHz maser internal Vlsr gradients appear to 
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Table 1. CH 3 OH 6.7 GHz EVN observations 



Source 


Epoch 1 


Epoch 2 


Epoch 3 


Beam FWHM 


Vel. Res. 


Ref. 










Maj. Axis 


Min. Axis 








(yr m d) 


(yr m d) 


(yr m d) 


(mas) 


(mas) 


(kms- 1 ) 




G16.59-0.05 


2006 Feb 26 


2007 Mar 16 


2008 Mar 15 


15 


5 


0.09 


SMC1 


G23.01-0.41 


2006 Feb 27 


2007 Mar 17 


2008 Mar 16 


13 


5 


0.09 


SMC2 


IRAS 20126+4104 


2004 Nov 06 


2007 Mar 21 


2009 Mar 11 


8 


6 


0.09 


MCR 


AFGL5142 


2004 Nov 04 


2007 Mar 16 


2009 Mar 12 


9 


6 


0.09 


GMS1, GMS2 



Notes. Column 1 gives the maser source name, and Columns 2, 3, and 4 the observing date of the first, second, and third EVN epochs, respectively. 
Columns 5 and 6 report the FWHM major and minor axes of the naturally-weighted beam used to reconstruct the maser images. Column 7 indicates 
the velocity resolution of the correlated visibilities, and Column 8 gives the bibliographic references. 



be remarkably stable in time on a timescale of (at least) several 
years. 

3.3. Ordered spatial and angular distribution of V L sr gradient 
directions 

For each of the four maser sources, Figure [5] shows the distribu- 
tion on the plane of the sky of the directions of Vlsr gradients 
and the histogram of the P.A. of gradient directions. To produce 
these plots, we neglected the sign of the Vlsr gradient, and the 
P.A. of negative gradients (180° < P s < 360°) was folded into 
the range 0° - 180°. 

Looking at the histograms of the gradient P.A., it is evident 
that, in each of the four sources, the directions of feature gradi- 
ents do not distribute uniformly, but tend to concentrate into spe- 
cific P.A. ranges. For the sources G23.01-0.41, G16.59-0.05, 
and IRAS 20126+4104, the P.A. distribution presents two broad, 
but clearly distinct peaks at the angles 40°+20° and 140°+20°. 
For the source AFGL 5142, only a single, broad peak, over the 
angles 90° +30°, is clearly visible. For the sources G23.0 1-0.41, 
G16.59-0.05, and IRAS 20126+4104, with a larger number of 
detected maser features and better statistics, the peaks in the dis- 
tribution of gradient P.A. become narrower when selecting gra- 
dients with increasingly better measured directions (i.e. corre- 
sponding to linear fits with increasingly higher values of the cor- 
relation coefficient, r s ). That reinforces the impression that the 
two-peak distribution has a physical meaning, and it does not 
result from a bias in our data analysis. 

The plots of the spatial distribution of the Vlsr gra- 
dient directions illustrate that, in the sources G23.0 1-0.41, 
G16.59-0.05, and IRAS 20126+4104, the two groups of fea- 
tures with gradients approximately perpendicular to each other, 
are not confined in two different regions of the sky, but are 
spread over a similar area and are observed to almost overlap 
(with separation of only ^10 mas) along some lines-of-sight. 
In these three maser sources, the two groups of gradients also 
have a similar spatial distribution of Vlsr- The case of the 
source AFGL 5142 may be different. In this source, the more 
red-shifted features, confined to the south and southeast of the 
whole maser distribution, have Vlsr gradients directed close to 
east-west (P.A. » 90°), whereas most of the "yellow" features, 
placed to the north and northwest, show gradients directed close 
to north-south (P.A. ss 0° or 180°). For this source, the poor 
statistics of well-measured gradient orientations very likely pre- 
vents us from clearly detecting a second peak in the distribu- 
tion of gradient P.A. close to 0° and/or 180°. In any case, in 
AFGL 5142, differently from the other three sources, the direc- 



tion of the maser Vlsr gradient appears to be correlated with the 
feature position. 

Section [5] discusses the peculiarities observed in the spatial 
and angular distribution of the 6.7 GHz maser Vlsr gradients, 
and proposes a simple kinematical interpretation. 

4. Comparison with 6.7 GHz maser proper motions 

The time persistency of most of the 6.7 GHz maser Vlsr gra- 
dients on timescales of (at least) several years, and the reg- 
ular (spatial and angular) distribution of the gradient direc- 
tions observed towards all the four maser targets suggests that 
the maser internal gradients might reflect the ordered gas mo- 
tions observed over the whole maser region (size of ~100- 
1000 AU) on small (linear) scales (~10 AU). In all sources but 
G16.59-0.05, our 22 GHz water maser observations show a fast 
wide-angle wind and/or a collimated jet emerging from the same 
(proto)star exciting the methanol masers. Over the whole region, 
the 6.7 GHz masers present a rather ordered 3-D kinematics, 
suggesting either rotation (in G16.59-0.05), a combination of 
rotation plus expansion (G23.01-0.41 and IRAS 20126+4104), 
orinfall (AFGL 5142). 

If the Vlsr gradients internal to the 6.7 GHz masers reflect 
ordered, large-scale gas motions, one would expect them to be 
related to the gas velocities. For instance, should the Vlsr gra- 
dients be produced in either a rotating structure or a Hubble 
outflow, one would expect the directions of the feature gradi- 
ent and proper motion to project at a close angle on the sky. 
Figure [6] presents the source-average distribution of the angle 
between the direction of the Vlsr gradient of a maser feature 
and its proper motion vector. The 6.7 GHz maser proper mo- 
tions are read from Table 3 of SMC1 for G16.59-0.05, from 
Table 4 of SMC2 for G23.0 1-0.41, from Table 1 of MCR for 
IRAS 20126+4104, and from GMS2 for AFGL 5142. The acute 
angle between the feature gradient and proper motion directions 
on the sky plane is used to measure the angular separation. To 
produce the plot in Fig. feature gradients with loosely defined 
directions (i.e. with r s < 0.5 and r v < 0.5) were discarded, and, 
for persistent features, the mean gradient averaged over the three 
epochs was calculated. Features that have both the Vlsr gradi- 
ent and the proper motion reliably measured, are a few tens for 
the sources G23.01-0.41 (36) and AFGL 5142 (20), but signif- 
icantly less for G16.59-0.05 (8) and IRAS 20126+4104 (7). To 
increase the statistics, we calculated the distribution of the frac- 
tional number of features for each source (dividing the number 
of features of a histogram bin by the total number of features), 
and then averaged over the four sources. 

Figure|6]provides some evidence that the Vlsr gradients tend 
to be oriented at relatively small angles (< 40°) from the direc- 
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tion of the maser proper motions. The actual distribution might 
be significantly more peaked at small angles than observed, since 
the uncertainty in the measurements of the gradient and proper 
motion orientations can be large (>20°), and can contribute sig- 
nificantly to broaden the peak. The basic characteristics of the 
source-average plot of Fig.[6]do not result from one source's dis- 
tribution dominating the other ones, but, conversely, in all four 
sources there is a similar tendency towards small angular separa- 
tions between the maser gradients and proper motions. Therefore 
we can conclude that in each of the observed sources, the data 
are consistent with having the Vlsr gradients and proper motion 
vectors in the same direction on the sky, considered the measure- 
me nt uncertainties. 

iFish et alJ ( 12005) . studying the 1.6 GHz OH masers towards 
the W3(OH) UC Hn region, did not note any correlation between 
the orientation of the Vlsr gr adients and the proper motions of 
the OH masers (measured by iBloemhof et all 1 1 9921). Over the 
whole exte nt of the UC Hn region (^300 AU in size). lFish et alJ 
(120061) and IFish & Siouwermanl (120071) did not find any large- 
scale organization in the distribution of the OH Vlsr gradients, 
and then suggest that the observed gradients should trace local 
phenomena, which are possibly associated to AU-scale turbu- 
lent fluctuations in the sky-projected velocity field. However, 
the environment traced by the 6.7 GHz masers in each of our 
target sources appears to be different from that of an UC Hn re- 
gion. The typical 6.7 GHz maser velocities (5-10 km s _1 ) are 
significantly higher than those of the OH masers in W3(OH) 
(a few km s _1 ). Whereas in W3(OH) and in other sources as 
well (G23. 01-0.41; see SMC2), the OH masers appear to mark 
gas at greater distance (> 10 3 -10 4 AU) from the (proto)star and 
to trace the slow expansion of the UC Hn region, the 6.7 GHz 
masers are more likely associated with gas accreting, rotating, 
and outflowing, in the proximity (within hundreds of AU) to the 
(proto)star. Differently from the OH, we believe that the mea- 
sured gradients of CH3OH masers trace gas kinematics on AU- 
scales, as discussed in detail in the next section. 

5. Nature of the 6.7 GHz maser V LS r gradients 

In the following, we discuss the hypothesis that the 6.7 GHz 
maser Vlsr gradients trace, on very small scale (<10 AU), the 
bulk-ordered motion measured on large scales from maser kine- 
matics: outflow or infall (Sect. 5.1) and rotation (Sect. 5.2) in 
high-mass (proto)stars. 

5.1. Do methanol maser V LSR gradients trace outflow/infall? 

We consider first the case that the 6.7 GHz maser Vlsr gradients 
are associated with the large-scale outflows traced by the wa- 
ter masers in the sources G23. 01-0.41, IRAS 20126+4104 and 
AFGL5142. 

In the source IRAS 20126+4104, the 3-D velocity field of the 
22 GHz water masers has been successfully reproduced with a 
model of a conical, Hubble flow, with maser velocities increas- 
ing linearly with the distance from the star. Assuming that the 
water maser velocity distribution can also be described in terms 
of an Hubble flow in the other sources, from the measured 3- 
D velocities it is possible to derive a typical velocity gradient in 
the outflowing gas. Considering that the observed water masers 
generally move close to the plane of the sky, one has to divide 
the maximum observed value of transverse velocities by the sky- 
projected semi-length of the outflow. This calculation gives a 
gradient of 0.2 km s _1 mas -1 and 0.15 km s _1 mas~' for the 
jets observed in the sources AFGL 5142 and IRAS 20126+4104, 



respectively. For the source G23. 01-0.41, we obtain a gradi- 
ent of 0.5 km s mas -1 for the fast and compact wide-angle 
wind located close to the (proto)star and a much lower value of 
0.06 km s _1 mas -1 for the collimated jet traced by two groups of 
water masers at a greater distance from the (proto)star. In each of 
the three sources, the gradient in the outflow velocities derived 
from the water masers is an upper limit to the average value of 
the 6.7 GHz maser Vlsr gradients (see Fig. [3]). That is consistent 
with the fact that 6.7 GHz masers are found to move much more 
slowly than water masers do. One possible interpretation is that 
water masers trace fast, shocked gas closer to the jet axis, while 
methanol masers might originate in relatively slow-moving ma- 
terial, entrained in the outflow at a larger distance from the jet 
axis. Looking at specific cases, in the source IRAS 20126+4104 
MCR show that the entraining process takes place for a group 
of 6.7 GHz features with proper motions directed close to the 
jet direction, while for the source AFGL 5142 we propose that 
the 6.7 GHz maser 3-D velocities trace infall rather than out- 
flow. Regardless of the specific source case, this simple discus- 
sion shows that, in each of the three sources where an outflow or 
infall is detected, the gas moves fast enough to be able to account 
for the Vlsr gradients internal to the 6.7 GHz masers. 

Figure [5] shows that the directions of the Vlsr gradients in 
three sources tend to cluster in two RA. ranges (separated by 
about 90°-100°), while a single peak in the distribution of gra- 
dient RA. is evident for the source AFGL 5142. For the three 
sources with an observed water maser jet, the orientation of the 
jet axis is reported in the plot of the spatial distribution of gra- 
dient directions. It is interesting to note that, for both sources 
G23.01-0.41 and IRAS 20126+4104, the PA. of the jet axis 
(about 50° and 115° for the first and second sources, respec- 
tively) is close to one of the two peaks in the distribution of gra- 
dient RA.. We consider that, while the direction of the jet axis is 
the one along which Vlsr gradients associated with the outflow 
should project, the direction perpendicular to the jet should be 
the preferred orientation of gradients due to (almost edge-on) ro- 
tation about the jet. This qualitative argument can simply explain 
the presence of the two peaks in the distribution of gradient RA. 
for the sources G23.01-0.41 and IRAS 20126+4104, and lends 
further support to the idea of a correlation between the small- 
scale Vlsr gradients inside the 6.7 GHz maser features and the 
ordered gas motions observed over the whole maser region. 

5.2. Do methanol maser V L sr gradients trace rotation? 

One of the simplest interpretations of a regular change in Vlsr 
with position across the plane of the sky is in terms of rotation. 
If the spots of a 6.7 GHz maser feature traced a portion of an an- 
nulus of a rotating sphere of gas, seen almost edge-on, that could 
explain the (sky-projected) linear distribution of spots and also 
the linear variation in Vlsr with (sky-projected) spot positions. 
Assuming that the rotation is gravitationally supported, one can 
easily derive a relation between the average gas (mass) density, 
p, and the measured gradient amplitude, F: 

r 2 =G^§ =4/3 nGp (1) 

where G is the gravitational constant, and M s and R s are 
the mass and the radius of the sphere, respectively. Using an 
average value of T = 0.05 km s _1 AU~' (as derived above 
for the sources IRAS 20126+4104 and G23.01-0.41), one finds 
p = 4 10~ 13 g crrT 3 , corresponding to a number density of hy- 
drogen molecules rtH 2 = 10 11 cm 4 . Using this average density, 
indicating with R\o the radius of the maser annulus in units of 
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10 AU, the mass of the gas sphere inside the annulus radius is 
M s = 2.7 l(T 3 /^ M . 

The obtained value of p and the expression for M s can be 
used to further constrain the properties of the putative rotation 
traced by the 6.7 GHz maser Vls r gradients. Excitation models 
of the 6.7 GHz methanol masers dCragg et al.ll2005l) predict that 
this maser emission should be quenched for n# 2 > 10 9 cm" 3 . 
Since this value is two orders of magnitude lower than the es- 
timate of the average density p of the masing sphere, we can 
conclude that it is unlikely that 6.7 GHz maser features trace a 
portion of an annulus of a rotating sphere of gas, unless the gas 
mass is very concentrated close to the center. Then, if the Vlsr 
gradients are tracing rotation, Keplerian rotation should be pre- 
ferred over solid-body rotation, which results from an homoge- 
nous mass distribution. The typical (sky-projected) size of a fea- 
ture's spot distribution (see examples reported in Fig.Q]) ranges 
from a few mas to 10 mas, or 5-50 AU, for a source distance 
from a few kpc up to 5 kpc. If the diameter of the maser an- 
nulus were comparable to the observed feature size, the mass 
M s enclosed inside the annulus radius, would be less or much 
less than a few percent of a solar mass. In each of the four 
sources, 6.7 GHz maser emission is observed close to a mas- 
sive proto(star) where the density and the temperature of the gas 
are relatively high {n Hl m 10 6 cm 3 , T* 100-200 K). Since the 
characteristic Jeans mass of such a dense and warm gas should 
be about 1 M , it is very unlikely that most of maser features in 
the observed sources trace a self-gravitating body with a mass 
of thousandths or hundredths of a solar mass. Therefore we fa- 
vor the alternative view where the annulus radius is significantly 
larger than the feature size, so that the enclosed mass M s would 
be comparable to or higher than 1 M . Our conclusion is that, 
if the gradients internal to the 6.7 GHz masers traced rotation, a 
maser feature (with size ~10 AU) would more likely trace only 
a small arc of a ring of matter in Keplerian rotation around a 
central mass ^ 1 M . 

In Section 15.2.11 we discuss the specific case of 
G16.59-0.05, where the methanol maser Vlsr gradients, under 
the assumption of Keplerian rotation, enabled us to constrain the 
(proto)star position and mass, finding good agreement with the 
values derived from the (3-D) maser spatial and velocity distri- 
bution. 



5.2.1. Vlsr gradients tracing rotation in G1 6.59-0.05 

Among the four maser targets, the source G16.59-0.05 is per- 
haps the one where the 3-D velocity field of the 6.7 GHz 
methanol masers traces the simplest kinematics, consisting in 
a rotating structure, elongated 073-074 towards the southeast- 
northwest direction and inclined ^30° with the plane of the 
sky (SMC1, Fig. 7). No collimated water maser jet is observed in 
this source, and it is likely that most of the 6.7 GHz Vlsr gradi- 
ents trace rotation, since the maser internal gradients are mainly 
directed transversally to the line connecting the maser position 
with the center of rotation (see Fig. [5J5). Now we apply the gen- 
eral arguments of Sect. I5.2I to this specific source, making use 
of the maser Vlsr gradients to constrain the (proto)star position 
and mass. Equation Q] indicates that, if maser internal gradients 
trace Keplerian rotation around a central mass M s , the product 
of the cube of the maser distance to the (proto)star, R 3 S , and the 
square of the gradient amplitude, F 2 , is a constant proportional 
to the stellar mass. Equation[T]is strictly valid for edge-on rota- 
tion. More in general, if the rotation axis is inclined by an angle 
i with the line-of-sight, for sufficiently low values of the polar 



angles, where < 9 < arctan(cos~ '(/)), the same relation holds 
replacing the central mass M s with the product M s sin 2 (/). The 
maser gradients in G16.59-0.05, tracing a rotating structure with 
an inclination angle i a:30°, should satisfy a relation equivalent 
to Eq.Q]over a wide range of polar angles (0 < 9 < 50°). 

Owing to the small inclination with the plane of the sky 
of the maser disk/toroid observed in G16.59-0.05, the sky- 
projected distances, R sp , should approximate the real 3-D dis- 
tances, R s , within 15%. Hence, we plotted the cube of the 
sky-projected maser distance to the (proto)star, R 3 sp , versus the 
square of the Vlsr gradient amplitude, F 2 , and fit to the data 
the curve R] p = K / F 2 , where K is a constant. In the as- 
sumption that the maser internal gradients trace Keplerian rota- 
tion, such a fit should produce minimal residuals when the sky- 
projected position of the (proto)star matches the real position. 
We searched for the protostar position over a grid centered on 
the 6.7 GHz maser center of motion (i.e. the (proto)star location 
independently estimated by averaging the position of the persis- 
tent maser features), with a grid semi-size of 200 mas at steps 
of 10 mas. Figure [7] presents the map of the fit residuals (mean 
squared) and best-fit for the (proto)stellar mass values. The posi- 
tion of the minimal fit residual ((East, Norfh)=(0,-0.02) arcsec) 
is found at an offset of only 20 mas (to the south) from the cen- 
ter of motion (at the origin of the maps shown in Fig. [7]). At the 
position of the minimal fit residual, the best-fit central mass (cor- 
rected for the inclination angle i = 30°) is 36 M , in optimal 
agreement with the value of 35 M derived from the 3-D veloc- 
ity field of the 6.7 GHz masers, assuming centrifugal equilibrium 
(SMC1). 

Thus, Fig.|7]shows that the maser internal gradients constrain 
the star position and mass to similar values as deduced from the 
maser spatial and velocity distributions. This result strongly sug- 
gests that, in G16.59-0.05, the 6.7 GHz maser gradients have a 
kinematical origin and are tracing the same (Keplerian) rotation 
marked by the large-scale maser kinematics. 

5.3. Caveats 

A more quantitative comparison between the ordered, large- 
scale motions traced with the (water and methanol) masers in 
each of the target sources, and the properties (direction and am- 
plitude) of the 6.7 GHz maser Vlsr gradients, is complicated by 
the limited knowledge of the maser geometry, and by the uncer- 
tainty in the position and mass of the (proto)star. Maser emis- 
sions often sample only a few directions, irregularly distributed 
across the (proto)stellar environment. Therefore, the (proto)star 
position on the plane of the sky, derived as the symmetry cen- 
ter of the maser spatial and velocity distribution, can be gener- 
ally determined with an accuracy not better than 50-100 mas 
(MCR). 

Furthermore, one cannot exclude that the maser Vlsr gra- 
dients reflect not just a single, large-scale ordered motion (as 
the simple cases of either outflow or rotation discussed in 
Sects. 15 . 1 1 and [572b , but a combination of them. That certainly 
would make it harder to interpret the maser internal gradients 
and to correlate them with the kinematics of the gas close to the 
(proto)star(s). The sources G23.01-0.41 and IRAS 20126+4104 
can be good examples of a rather complex kinematical pattern 
of the 6.7 GHz methanol masers, which trace rotation (about 



2 The position of the orbiting maser feature is described using polar 
coordinates centered on the (proto)star, with the polar angle 6 taken to 
be zero when the polar radius is at the minimal separation (90 - z) from 
the line-of-sight. 
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the protostellar jet) and outflow (across and along the jet). 
Assessment of the potential of the new tool provided by the 
6.7 GHz V L sr gradients to constrain the gas kinematics would 
certainly benefit from an increased statistics of well-studied 
high-mass star-forming regions. 

6. Conclusion 

This work presents a study at high-velocity resolution of the mil- 
liarcsecond structure of the 6.7 GHz CH3OH masers observed 
towards four high-mass star-forming regions: G16.59-0.05, 
G23.0 1-0.41, IRAS 20126+4104, and AFGL 5142. Most of the 
detected 6.7 GHz maser features present an ordered (linear or 
arc-like) distribution of maser spots on the plane of the sky, to- 
gether with a regular variation in the spot Vlsr with position 
measured along the major axis of the feature elongation. A fea- 
ture Vlsr gradient can be defined as a vector quantity, character- 
ized by an orientation (that of the feature major axis) and by an 
amplitude (equal to the derivative of the spot V L sr with position). 
Typical values for the amplitude of the 6.7 GHz maser Vlsr gra- 
dients are found to be 0.1-0.2 km s _1 mas -1 . Our multi-epoch 
VLBI observations show that, in each of the four target sources, 
the orientation and the amplitude of most of the feature Vlsr 
gradients remain remarkably stable in time, on timescales of (at 
least) several years. In three (G16.59-0.05, G23.01-0.41, and 
IRAS 20126+4104) of the four sources under examination, the 
directions of the feature gradients on the plane of the sky con- 
centrate in two intervals of R A., each s;30 o ^10 o wide, separated 
by about 90°-100°. The two groups of features with V L sr gradi- 
ents oriented approximately perpendicular to each other, present 
a similar spatial distribution over the whole maser region. 

The observed time persistency of the 6.7 GHz maser V L sr 
gradients and the order found in their angular and spatial dis- 
tributions suggests a kinematical interpretation for their origin. 
That is supported further by the finding that, in each of the four 
sources, the data are consistent with having the V L sr gradients 
and proper motion vectors in the same direction on the sky, con- 
sidered the measurement uncertainties. We discuss the case that 
the observed V L sr gradients can reflect the ordered, large-scale 
(~100-1000 AU) motions (outflow, rotation, and infall) traced 
by the (22 GHz water and 6.7 GHz methanol) masers in these 
regions, on much smaller linear (~10 AU) scales. For the source 
G16.59-0.05, where the 6.7 GHz masers trace a well-defined ro- 
tating structure, we have used the maser V L sr gradients to con- 
strain the (proto)star position and mass and find good agreement 
with the values derived from the (3-D) maser spatial and velocity 
distribution. 

Our conclusion is that the study of the 6.7 GHz maser gradi- 
ents on milliarcsecond scales proves to be a useful tool for inves- 
tigating the gas kinematics in proximity to massive (proto)stars. 
In the future, we plan to extend the analysis of the 6.7 GHz maser 
Vlsr gradients to other sources, with the particular aim of fur- 
ther investigating their origin and, ultimately, employing them 
for kinematical studies. 
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(a) 



km s- 1 




Major Axis offset (mas) 

Fig. 1. Examples of two 6.7 GHz maser features with ordered spatial and Vlsr distribution, (a) The upper, middle, and lower panels refer to measurements at the 
first, second, and third observing epochs, respectively. Left panels: Spot spatial distribution. Triangles, squares, and pentagons indicate the spot positions on the 
plane of the sky at the first, second, and third epochs, respectively, with symbol size varying logarithmically with the spot intensity. The symbol color indicates the 
spot Vlsr, with the color-velocity conversion code given in the wedge on the right side of the three panels. At each epoch, the positions of the spots are relative to 
that of the epoch's most intense spot. In each of the three panels, the continuous line shows the (minimum^ 2 ) linear fit to the spot positions, defining the direction of 
the major axis of the spot distribution. The lower left corner of each panel reports the value of the linear correlation coefficient, r s . Right panels: Spot VLSR-position 
distribution. In each of the three panels, empty circles show the plot of spot Vlsr versus the position offset along the major axis of the spot distribution. Position 
offsets are taken as positive if increasing towards the east. Circle size varies logarithmically with the spot intensity. The continuous lines show the (minimum^ 2 ) 
linear fit to the circle distribution. The lower right corner of each panel reports the value of the linear correlation coefficient, r v . (b) Same as for (a) for another 
feature of our sample. 
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Fig. 2. Histograms of the linear correlation coefficients r s and r v . Each panel presents the histograms of r s (blue dashed line) and 
r v (red dashed line), for one of the four maser sources, indicated in the upper left corner of the panel. The bin size is 0. 1 . 
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Fig. 3. Histograms of the amplitude of the V L sr gradient, F. Each panel plots data of a single maser source, indicated in the upper 
right corner of the panel. The bin size of all the histograms is 0.05 km s _1 mas 1 . Each panel (below the source name) also reports 
the arithmetic average, fi, and the standard deviation, cr, of the gradient amplitude, in units of [km s~' mas -1 ]. 
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Fig. 4. Time variation in Vlsr gradient P.A. and amplitude, (a) Each of the four panels refers to the maser source indicated in the upper right 
corner of the panel, and presents the histogram of the standard deviation (A P s ) of the values of the gradient P.A., at different epochs, for features 
persisting over two or three observing epochs. The bin size of all the histograms is 10°. (b) Each of the four panels refers to the maser source 
indicated in the upper right corner of the panel, and presents the histogram of the fractional time variation (A T / T) in the gradient amplitude, for 
features persisting over two or three observing epochs. The bin size of all the histograms is 0.2. 
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(b) G16.59-0.05 
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Fig. 5. The panels in each row present the spatial distribution and the histogram of the Vlsr gradient directions for the maser source indicated on the left side above 
the panels. Left Panel: Spatial distribution of the Vlsr gradient directions. Colored dots report the position of the maser features detected over the three observing 
epochs, with colors indicating the feature Vlsr- The velocity-color conversion code is shown in the wedge on the right side of the panel, with green denoting the 
systemic Vlsr of the maser source. Colored segments associated to maser features give the direction of the feature Vlsr gradient (with P.A. varying in the range 
0°-180°), with different colors to distinguish the observing epoch: red, green, and blue, for the first, second, and third epochs, respectively. The segment length is 
proportional to the value of the correlation coefficient,^, of the linear fit to the spot positions, from which the gradient direction is derived. The plot reports only 
the feature gradients with better defined directions, corresponding to more linear spot distributions and higher values of the correlation coefficient, r s : r s > 0.9 for 
G23.01-0.41, and r s > 0.8 for G16.59-0.05. For the sources G16.59-0.05 and G23.01— 0.41, feature positions are relative to the maser's "center of motion", as 
defined in SMCI and SMC2, respectively. In the plot of G23.01-0.41, the dashed arrow indicates the direction of the collimated jet traced close to the (proto)star 
by the 22 GHz water masers (SMC2). Right panel: Histograms of the P.A. of the feature gradient directions. Dotted black, dashed red, and continuous blue lines 
show the histograms of the gradient P.A. for features with increasingly better linear structure, corresponding to values of the correlation coefficient r s > Q,r s > 0.5, 
and r s > 0.9, respectively. The histogram bin size is 15° for the source G23. 01-0.41, and 20° for the source G16.59-0.05. 
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Fig. 5. The panels in each row present the spatial distribution and the histogram of the Vlsr gradient directions for the maser source indicated on the left side above 
the panels. Left panel: Spatial distribution of the Vlsr gradient directions. Colored dots and colored segments have the same meaning as in Figs.fS^ andfSp. For the 
source IRAS 20126+4104, absolute maser positions are reported, relative to the observing epoch 2004 November 6 (M CR, Fig. 2). In the source AFGL 5142, feature 
positions are relative to the peak of the VLA 1.3 cm continuum observed toward the 6.7 GHz masers (Goddi et al. 2007, see Fig. 2), whose contour map is shown by 
a dotted line. In both the plots of IRAS 20126+4104 and AFGL 5142, the dashed line gives the direction of the collimated jet traced close to the (proto)star by the 
22 GHz water masers. The plot reports only the feature gradients with better defined directions, corresponding to more linear spot distributions and higher values of 
the correlation coefficient, r s : r s > 0.8 for IRAS 20126+4104, and r s > 0.6 for AFGL 5 142. Right panel: Histograms of the RA. of the feature gradient directions. 
Dotted black, dashed red, and continuous blue lines have the same meaning as in Figs. [5^ and[5p. The histogram bin size is 15° for the source IRAS 20126+4104, 
and 20° for the source AFGL 5142. 
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Fig. 6. Histogram of the distribution of the angle between the direction of the Vlsr gradient of a maser feature and its proper motion 
vector. The histogram bin size is 20°. The histogram reports the fractional number of maser features averaged over the four target 
sources. 





Fig. 7. Maps of the (mean squared) residual (left panel) and best-fit central mass (right panel) obtained by fitting the curve 
R] p — K I T 2 to the 6.7 GHz maser data in G16.59-0.05. Colours give the map value, with the value-color conversion code given in 

the wedge on the right side of the panels. The values reported in the wedge are in units of arcsec 6 and solar mass (M ) for the left 
and right panels, respectively. In both panels, positions are given relative to the center of motion of the 6.7 GHz masers (as defined 
in SMC1). To better compare the two maps, in both panels white contours report the same levels of fit (mean squared) residual. The 
plotted levels are from 7 1CT 6 arcsec 6 (10% higher than the minimum fit residual) to 1.6 10~ 3 arcsec 6 , stepped by a factor of 2. The 
spatial distribution of the 6.7 GHz maser features is indicated with white dots. 
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